Using the enrichment analysis to compare the rankings and metricsf
#install.packages("rlist")
library(rlist)
Attaching package: ‘rlist’
The following object is masked from ‘package:S4Vectors’:
List
library(ggplot2)
library(latex2exp)
library(systemPipeR)
setwd("~/Dropbox/PNNL/Omics/SARS-COV-2/SARSCov2/")
source("GeneSets.R")
# get the background gene set/list
biggerTransGeneSet <- as.character(read.table("biggerTransGenes.txt", header = FALSE, sep = "")$V1)
bigTransGeneSet <- as.character(read.table("bigTransGenes.txt", header = FALSE, sep = "")$V1)
bakgrdGeneList <- bigTransGeneSet
lenBakgrd <- length(bakgrdGeneList)
#
# sars1Targets <- as.character(read.table("sars1-targets.txt", header = FALSE, sep = "")$V1)
# sars2Targets <- as.character(read.table("sars1-targets.txt", header = FALSE, sep = "")$V1)
# sarsTargets <- as.character(read.table("sars-targets.txt", header = FALSE, sep = "")$V1)
esFn <- function(x, p = 1) {
len = length(x)
miss = ifelse(x!=0, 0, 1)
#hiti = ((miss * 2 - 1) - miss) / 2
if (p == 0) {
xp <- ifelse(x!=0, 1, 0)
} else {
xp = abs(x)**p
}
sumxp = sum(xp)
summiss = sum(miss)
if (sumxp == 0) {
hitp = rep(0, len)
} else {
hitp = cumsum(xp) / sumxp
}
if (summiss == 0) {
missp = rep(0, len)
} else {
missp = cumsum(miss) / summiss
}
maximum = max(hitp - missp)
minimum = min(hitp - missp)
score = ifelse(abs(maximum) > abs(minimum), maximum, minimum)
return(score)
}
library(MASS)
epFn <- function(x, xnull) {
fit <- fitdistr(xnull[xnull > 0], densfun="normal")
pvalue = 1 - pnorm(abs((x - fit$estimate[1]) / (fit$estimate[2])))
return(pvalue)
}
Use the list of gene to do enrichment analysis
samplesize = 1000
targetSetNames <- c("sars1", "sars", "influenza", "influenzaIntersect", "influenzaAll", "ebola", "ebolaLong")
centralities <- c("s-closeness", "s-betweenness")
rankMethods <- c("Number of Increases", "Euclidean Distance", "Absolute Rank Change", "Relative Rank Change", "Absolute Change", "Relative Change", "Combined Score")
# rankMethods <- c("Number.of.Increases", "Euclidean.Distance", "Absolute.Rank.Change", "Relative.Rank.Change", "Absolute.Change", "Relative.Change", "Combined.Score")
virusNames <- c("Sars", "Mers", "Iv", "Eb")
```r
dfEnrich <- setNames(data.frame(matrix(ncol = 8, nrow = 0)), c(\genecount\, \p\, \Enrichment score\, \p-value\, \Rank method\, \Virus name\, \Target sets\,\Hypergraph metrics\))
dfi <- 1
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
### analyze the SARS1 targets with all ranks and metrics
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuZm9yIChjZW50cmFsaXR5IGluIGNlbnRyYWxpdGllcykge1xuICBmb3IgKHZpcnVzTmFtZSBpbiB2aXJ1c05hbWVzKSB7XG4gICAgaW5maWxlIDwtIHNwcmludGYoXFwuL25vJXMtJXNfY2hhbmdlcy5jc3ZcXCwgdmlydXNOYW1lLCBjZW50cmFsaXR5KVxuICAgIGlmIChmaWxlLmV4aXN0cyhpbmZpbGUpKSB7XG4gICAgICBnZW5lbGlzdCA8LSByZWFkLmNzdihpbmZpbGUsIGhlYWRlciA9IFRSVUUsIHJvdy5uYW1lcyA9IDEpXG4gICAgICBjb2xuYW1lcyhnZW5lbGlzdCkgPC0gcmFua01ldGhvZHNcbiAgICAgIGFsbEdlbmVuYW1lcyA8LSByb3duYW1lcyhnZW5lbGlzdClcbiAgICAgIGxlbkFsbCA8LSBsZW5ndGgoYWxsR2VuZW5hbWVzKVxuICAgICAgXG4gICAgICBmb3IgKHJhbmtNZXRob2QgaW4gcmFua01ldGhvZHMpIHtcbiAgICAgICAgcmFua09yaSA8LSBhcy5tYXRyaXgoZ2VuZWxpc3Rbb3JkZXIoLWdlbmVsaXN0W3JhbmtNZXRob2RdKSxdW3JhbmtNZXRob2RdKVxuICAgICAgICByYW5rIDwtIHN1YnNldChyYW5rT3JpLCBpcy5maW5pdGUocmFua09yaSkpXG4gICAgICAgIGdlbmVuYW1lcyA8LSByb3duYW1lcyhyYW5rKVxuICAgICAgICBcbiAgICAgICAgZm9yICh0YXJnZXRTZXROYW1lIGluIHRhcmdldFNldE5hbWVzKSB7XG4gICAgICAgICAgdGFyZ2V0U2V0IDwtIGFzLmNoYXJhY3RlcihyZWFkLnRhYmxlKHBhc3RlMCh0YXJnZXRTZXROYW1lLFxcLXRhcmdldHMudHh0XFwpLCBoZWFkZXIgPSBGQUxTRSwgc2VwID0gXFxcXCkkVjEpXG4gICAgICAgICAgZW5yIDwtIGVucmljaG1lbnRfYnlfZmlzaGVycyhnZW5lbmFtZXMsIGJha2dyZEdlbmVMaXN0LCB0YXJnZXRTZXQpXG4gICAgICAgICAgcCA9IGVuciRmaXNoZXIkcC52YWx1ZVxuICAgICAgICAgIGYgPSBlbnIkZm9sZHhcbiAgICAgICAgICBtYXQgPSBlbnIkbWF0XG4gICAgICAgICAgbGVuU2hvcnQgPC0gbGVuZ3RoKGdlbmVuYW1lcylcbiAgICAgICAgICB6c2NvcmVzIDwtIHJhbmtcbiAgICAgICAgICB6c2NvcmVzW3doaWNoKCEoZ2VuZW5hbWVzICVpbiUgdGFyZ2V0U2V0KSldIDwtIDBcbiAgICAgICAgICBlcyA8LSBlc0ZuKHpzY29yZXMpXG4gICAgICAgICAgc2FtcGxlZFZhbHVlcyA9IGMoKVxuICAgICAgICAgIGZvciAoayBpbiAxOnNhbXBsZXNpemUpIHtcbiAgICAgICAgICAgIHNhbXBsZWRpcyA9IHNhbXBsZSgxOmxlbkJha2dyZClcbiAgICAgICAgICAgIHRlbXBHZW5lbmFtZXMgPC0gZ2VuZW5hbWVzW3NhbXBsZWRpc1sxOmxlblNob3J0XV1cbiAgICAgICAgICAgIHRlbXBac2NvcmVzIDwtIHpzY29yZXNbc2FtcGxlZGlzWzE6bGVuU2hvcnRdXVxuICAgICAgICAgICAgdGVtcFpzY29yZXNbd2hpY2goISh0ZW1wR2VuZW5hbWVzICVpbiUgdGFyZ2V0U2V0KSldIDwtIDBcbiAgICAgICAgICAgIHNhbXBsZWRWYWx1ZXMgPSBjKHNhbXBsZWRWYWx1ZXMsIGVzRm4odGVtcFpzY29yZXMpKVxuICAgICAgICAgIH1cbiAgICAgICAgICBlcCA9IChsZW5ndGgod2hpY2goc2FtcGxlZFZhbHVlcyA+PSBlcykpICsgMSkgLyAoc2FtcGxlc2l6ZSArIDEpXG4gICAgICAgICAgI2VwID0gZXBGbihlcywgc2FtcGxlZFZhbHVlcylcbiAgICAgICAgICAjIHNhdmUgdGVzdCByZXN1bHRcbiAgICAgICAgICBkZkVucmljaFtkZmksXSA8LSBjKG1hdFsxLDFdLCBwLCBlcywgZXAsIHJhbmtNZXRob2QsIHZpcnVzTmFtZSwgdGFyZ2V0U2V0TmFtZSwgY2VudHJhbGl0eSlcbiAgICAgICAgICBkZmkgPC0gZGZpICsgMVxuICAgICAgICB9XG4gICAgICB9XG4gICAgfVxuICB9XG59XG5gYGBcbmBgYCJ9 -->
```r
```r
for (centrality in centralities) {
for (virusName in virusNames) {
infile <- sprintf(\./no%s-%s_changes.csv\, virusName, centrality)
if (file.exists(infile)) {
genelist <- read.csv(infile, header = TRUE, row.names = 1)
colnames(genelist) <- rankMethods
allGenenames <- rownames(genelist)
lenAll <- length(allGenenames)
for (rankMethod in rankMethods) {
rankOri <- as.matrix(genelist[order(-genelist[rankMethod]),][rankMethod])
rank <- subset(rankOri, is.finite(rankOri))
genenames <- rownames(rank)
for (targetSetName in targetSetNames) {
targetSet <- as.character(read.table(paste0(targetSetName,\-targets.txt\), header = FALSE, sep = \\)$V1)
enr <- enrichment_by_fishers(genenames, bakgrdGeneList, targetSet)
p = enr$fisher$p.value
f = enr$foldx
mat = enr$mat
lenShort <- length(genenames)
zscores <- rank
zscores[which(!(genenames %in% targetSet))] <- 0
es <- esFn(zscores)
sampledValues = c()
for (k in 1:samplesize) {
sampledis = sample(1:lenBakgrd)
tempGenenames <- genenames[sampledis[1:lenShort]]
tempZscores <- zscores[sampledis[1:lenShort]]
tempZscores[which(!(tempGenenames %in% targetSet))] <- 0
sampledValues = c(sampledValues, esFn(tempZscores))
}
ep = (length(which(sampledValues >= es)) + 1) / (samplesize + 1)
#ep = epFn(es, sampledValues)
# save test result
dfEnrich[dfi,] <- c(mat[1,1], p, es, ep, rankMethod, virusName, targetSetName, centrality)
dfi <- dfi + 1
}
}
}
}
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuZm9yIChjZW50cmFsaXR5IGluIGNlbnRyYWxpdGllcykge1xuICByYW5rTWV0aG9kIDwtIFxcaW5kaXZpZHVhbFxcXG4gIGluZmlsZSA8LSBzcHJpbnRmKFxcLi9iaWdUcmFucy0lcy0lcy1yYW5rcy5jc3ZcXCwgY2VudHJhbGl0eSwgcmFua01ldGhvZClcbiAgaWYgKGZpbGUuZXhpc3RzKGluZmlsZSkpIHtcbiAgICBnZW5lbGlzdCA8LSByZWFkLmNzdihpbmZpbGUsIGhlYWRlciA9IFRSVUUsIHJvdy5uYW1lcyA9IDEpXG4gICAgYWxsR2VuZW5hbWVzIDwtIHJvd25hbWVzKGdlbmVsaXN0KVxuICAgIGxlbkFsbCA8LSBsZW5ndGgoYWxsR2VuZW5hbWVzKVxuICAgIGZvciAodmlydXNOYW1lIGluIHZpcnVzTmFtZXMpIHtcbiAgICAgIHJhbmtPcmkgPC0gYXMubWF0cml4KGdlbmVsaXN0W29yZGVyKC1nZW5lbGlzdFt2aXJ1c05hbWVdKSxdW3ZpcnVzTmFtZV0pXG4gICAgICByYW5rIDwtIHN1YnNldChyYW5rT3JpLCBpcy5maW5pdGUocmFua09yaSkpXG4gICAgICBnZW5lbmFtZXMgPC0gcm93bmFtZXMocmFuaylcbiAgICAgIGZvciAodGFyZ2V0U2V0TmFtZSBpbiB0YXJnZXRTZXROYW1lcykge1xuICAgICAgICB0YXJnZXRTZXQgPC0gYXMuY2hhcmFjdGVyKHJlYWQudGFibGUocGFzdGUwKHRhcmdldFNldE5hbWUsXFwtdGFyZ2V0cy50eHRcXCksIGhlYWRlciA9IEZBTFNFLCBzZXAgPSBcXFxcKSRWMSlcbiAgICAgICAgZW5yIDwtIGVucmljaG1lbnRfYnlfZmlzaGVycyhnZW5lbmFtZXMsIGJha2dyZEdlbmVMaXN0LCB0YXJnZXRTZXQpXG4gICAgICAgIHAgPSBlbnIkZmlzaGVyJHAudmFsdWVcbiAgICAgICAgZiA9IGVuciRmb2xkeFxuICAgICAgICBtYXQgPSBlbnIkbWF0XG4gICAgICAgIGxlblNob3J0IDwtIGxlbmd0aChnZW5lbmFtZXMpXG4gICAgICAgIHpzY29yZXMgPC0gcmFua1xuICAgICAgICB6c2NvcmVzW3doaWNoKCEoZ2VuZW5hbWVzICVpbiUgdGFyZ2V0U2V0KSldIDwtIDBcbiAgICAgICAgZXMgPC0gZXNGbih6c2NvcmVzKVxuICAgICAgICBzYW1wbGVkVmFsdWVzID0gYygpXG4gICAgICAgIGZvciAoayBpbiAxOnNhbXBsZXNpemUpIHtcbiAgICAgICAgICBzYW1wbGVkaXMgPSBzYW1wbGUoMTpsZW5CYWtncmQpXG4gICAgICAgICAgdGVtcEdlbmVuYW1lcyA8LSBnZW5lbmFtZXNbc2FtcGxlZGlzWzE6bGVuU2hvcnRdXVxuICAgICAgICAgIHRlbXBac2NvcmVzIDwtIHpzY29yZXNbc2FtcGxlZGlzWzE6bGVuU2hvcnRdXVxuICAgICAgICAgIHRlbXBac2NvcmVzW3doaWNoKCEodGVtcEdlbmVuYW1lcyAlaW4lIHRhcmdldFNldCkpXSA8LSAwXG4gICAgICAgICAgc2FtcGxlZFZhbHVlcyA9IGMoc2FtcGxlZFZhbHVlcywgZXNGbih0ZW1wWnNjb3JlcykpXG4gICAgICAgIH1cbiAgICAgICAgZXAgPSAobGVuZ3RoKHdoaWNoKHNhbXBsZWRWYWx1ZXMgPj0gZXMpKSArIDEpIC8gKHNhbXBsZXNpemUgKyAxKVxuICAgICAgICAjZXAgPSBlcEZuKGVzLCBzYW1wbGVkVmFsdWVzKVxuICAgICAgICAjIHNhdmUgdGVzdCByZXN1bHRcbiAgICAgICAgZGZFbnJpY2hbZGZpLF0gPC0gYyhtYXRbMSwxXSwgcCwgZXMsIGVwLCByYW5rTWV0aG9kLCB2aXJ1c05hbWUsIHRhcmdldFNldE5hbWUsIGNlbnRyYWxpdHkpXG4gICAgICAgIGRmaSA8LSBkZmkgKyAxXG4gICAgICB9XG4gICAgfVxuICB9XG59XG5cbmBgYFxuYGBgIn0= -->
```r
```r
for (centrality in centralities) {
rankMethod <- \individual\
infile <- sprintf(\./bigTrans-%s-%s-ranks.csv\, centrality, rankMethod)
if (file.exists(infile)) {
genelist <- read.csv(infile, header = TRUE, row.names = 1)
allGenenames <- rownames(genelist)
lenAll <- length(allGenenames)
for (virusName in virusNames) {
rankOri <- as.matrix(genelist[order(-genelist[virusName]),][virusName])
rank <- subset(rankOri, is.finite(rankOri))
genenames <- rownames(rank)
for (targetSetName in targetSetNames) {
targetSet <- as.character(read.table(paste0(targetSetName,\-targets.txt\), header = FALSE, sep = \\)$V1)
enr <- enrichment_by_fishers(genenames, bakgrdGeneList, targetSet)
p = enr$fisher$p.value
f = enr$foldx
mat = enr$mat
lenShort <- length(genenames)
zscores <- rank
zscores[which(!(genenames %in% targetSet))] <- 0
es <- esFn(zscores)
sampledValues = c()
for (k in 1:samplesize) {
sampledis = sample(1:lenBakgrd)
tempGenenames <- genenames[sampledis[1:lenShort]]
tempZscores <- zscores[sampledis[1:lenShort]]
tempZscores[which(!(tempGenenames %in% targetSet))] <- 0
sampledValues = c(sampledValues, esFn(tempZscores))
}
ep = (length(which(sampledValues >= es)) + 1) / (samplesize + 1)
#ep = epFn(es, sampledValues)
# save test result
dfEnrich[dfi,] <- c(mat[1,1], p, es, ep, rankMethod, virusName, targetSetName, centrality)
dfi <- dfi + 1
}
}
}
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuc2F2ZShkZkVucmljaCwgZmlsZSA9IFxcZGZFbnJpY2guUkRhdGFcXClcbmBgYFxuYGBgIn0= -->
```r
```r
save(dfEnrich, file = \dfEnrich.RData\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxubG9hZChcImRmRW5yaWNoLlJEYXRhXCIpXG5gYGAifQ== -->
```r
load("dfEnrich.RData")
dfEnrich$`Enrichment score` <- as.numeric(dfEnrich$`Enrichment score`)
dfEnrich$`p-value` <- as.numeric(dfEnrich$`p-value`)
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "sars1",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "influenza",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "influenzaRNAi1",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "influenzaAll",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "ebola",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
# df <- dfEnrich[Reduce(intersect, list(which(dfEnrich$`Target sets` == "sars1"), which(), which(), which())),]
df <- dfEnrich[dfEnrich$`Target sets` == "ebolaLong",]
ggplot(data = df, mapping = aes(x=`Rank method`, y=`Enrichment score`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_continuous() + theme(axis.text.x = element_text(angle=270))
ggplot(data = df, mapping = aes(x=`Rank method`, y=`p-value`, color = `Virus name`, shape = `Hypergraph metrics`)) + geom_point(size=3 ,alpha = 0.6) + geom_line(linetype = "dashed") + scale_y_log10() + theme(axis.text.x = element_text(angle=270))
```r
targetSetNames0 <- targetSetNames <- c(\sars1\, \sars2\, \influenza\, \ebola\)
targetSets <- list()
for (i in seq(length(targetSetNames0))) {
targetSet <- as.character(read.table(paste0(targetSetNames0[i],\-targets.txt\), header = FALSE, sep = \\)$V1)
targetSets[[i]] <- targetSet
}
names(targetSets) <- targetSetNames0
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxudmVubnNldCA8LSBvdmVyTGFwcGVyKHRhcmdldFNldHMsIHR5cGU9XFx2ZW5uc2V0c1xcKVxudmVublBsb3QodmVubnNldCwgbXltYWluID0gXFxJbnRlcnNlY3Rpb25zIG9mIHRhcmdldCBnZW5lc1xcKVxuYGBgXG5gYGAifQ== -->
```r
```r
vennset <- overLapper(targetSets, type=\vennsets\)
vennPlot(vennset, mymain = \Intersections of target genes\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxudmVubkZpbGUgPC0gcGFzdGUwKFxcdmVubkRpYWdyYW0ucGRmXFwpXG5wZGYodmVubkZpbGUsaGVpZ2h0PTYsd2lkdGg9OClcbnZlbm5QbG90KHZlbm5zZXQsIG15bWFpbiA9IFxcSW50ZXJzZWN0aW9ucyBvZiB0YXJnZXQgZ2VuZXNcXClcbmRldi5vZmYoKVxuYGBgXG5gYGAifQ== -->
```r
```r
vennFile <- paste0(\vennDiagram.pdf\)
pdf(vennFile,height=6,width=8)
vennPlot(vennset, mymain = \Intersections of target genes\)
dev.off()
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuaW50ZXJMaXN0IDwtIHZlbm5saXN0KHZlbm5zZXQpXG5zYXZlKGludGVyTGlzdCwgZmlsZSA9IFxcaW50ZXJzZWN0aW9ucy5SRGF0YVxcKVxuYGBgXG5gYGAifQ== -->
```r
```r
interList <- vennlist(vennset)
save(interList, file = \intersections.RData\)
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxuaW50ZXJOYW1lcyA8LSBuYW1lcyhpbnRlckxpc3QpXG5mb3IgKGkgaW4gc2VxKGxlbmd0aChpbnRlck5hbWVzKSkpIHtcbiAgaW50ZXJOYW1lIDwtIGludGVyTmFtZXNbaV1cbiAgaW50ZXJnZW5lcyA8LSBpbnRlckxpc3RbW2ldXVxuICBpZiAobGVuZ3RoKGludGVyZ2VuZXMpID4gMCkge1xuICAgIHdyaXRlLnRhYmxlKGludGVyZ2VuZXMsIGZpbGUgPSBwYXN0ZTAoXFx0YXJnZXRTZXRzL2ludGVyc2VjdGlvbnMtXFwsIGludGVyTmFtZSwgXFwudHh0XFwpLCBzZXA9XFxcXHRcXCwgcXVvdGU9RkFMU0UsIGNvbC5uYW1lcyA9IEZBTFNFKVxuICB9XG59XG5gYGBcbmBgYCJ9 -->
```r
```r
interNames <- names(interList)
for (i in seq(length(interNames))) {
interName <- interNames[i]
intergenes <- interList[[i]]
if (length(intergenes) > 0) {
write.table(intergenes, file = paste0(\targetSets/intersections-\, interName, \.txt\), sep=\\t\, quote=FALSE, col.names = FALSE)
}
}
<!-- rnb-source-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxudG9wTnVtIDwtIDEwMDBcbmBgYCJ9 -->
```r
topNum <- 1000
```r
for (rankMethod in rankMethods) {
for (centrality in centralities) {
targetSets <- list()
for (i in seq(length(virusNames))) {
infile <- sprintf(\./no%s-%s_changes.csv\, virusNames[i], centrality)
if (file.exists(infile)) {
genelist <- read.csv(infile, header = TRUE, row.names = 1)
colnames(genelist) <- rankMethods
rankOri <- as.matrix(genelist[order(-genelist[rankMethod]),][rankMethod])
rank <- subset(rankOri, is.finite(rankOri))
genenames <- rownames(rank)
targetSets[[i]] <- genenames[1:topNum]
}
}
names(targetSets) <- virusNames
label = paste0(centrality, \-\, rankMethod, \-top-\, topNum)
dir.create(file.path(\venn\, label))
vennset <- overLapper(targetSets, type=\vennsets\)
vennFile <- paste0(\venn/\, label, \.pdf\)
pdf(vennFile,height=6,width=8)
vennPlot(vennset, mymain = paste0(\Intersections of target genes: \, label))
dev.off()
interList <- vennlist(vennset)
save(interList, file = paste0(\venn/\, label, \.RData\))
interNames <- names(interList)
for (i in seq(length(interNames))) {
interName <- interNames[i]
intergenes <- interList[[i]]
if (length(intergenes) > 0) {
write.table(intergenes, file = paste0(\venn/\, label, \/\, interName, \.txt\), sep=\\t\, quote=FALSE, col.names = FALSE, row.names = FALSE)
}
}
}
}
<!-- rnb-source-end -->
<!-- rnb-output-begin eyJkYXRhIjoiJ3Zlbm4vcy1jbG9zZW5lc3MtTnVtYmVyIG9mIEluY3JlYXNlcy10b3AtMTAwMCcgYWxyZWFkeSBleGlzdHNcbiJ9 -->
‘venn/s-closeness-Number of Increases-top-1000’ already exists
<!-- rnb-output-end -->
<!-- rnb-chunk-end -->
<!-- rnb-text-begin -->
<!-- rnb-text-end -->
<!-- rnb-chunk-begin -->
<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuYGBgclxucmFua01ldGhvZCA8LSBcXGluZGl2aWR1YWxcXFxuXG5mb3IgKGNlbnRyYWxpdHkgaW4gY2VudHJhbGl0aWVzKSB7XG4gIFxuICB0YXJnZXRTZXRzIDwtIGxpc3QoKVxuICBpbmZpbGUgPC0gc3ByaW50ZihcXC4vYmlnVHJhbnMtJXMtJXMtcmFua3MuY3N2XFwsIGNlbnRyYWxpdHksIHJhbmtNZXRob2QpXG4gIGZvciAoaSBpbiBzZXEobGVuZ3RoKHZpcnVzTmFtZXMpKSkge1xuXG4gICAgaWYgKGZpbGUuZXhpc3RzKGluZmlsZSkpIHtcbiAgICAgIGdlbmVsaXN0IDwtIHJlYWQuY3N2KGluZmlsZSwgaGVhZGVyID0gVFJVRSwgcm93Lm5hbWVzID0gMSlcbiAgICAgIGNvbG5hbWVzKGdlbmVsaXN0KSA8LSB2aXJ1c05hbWVzXG4gICAgICBcbiAgICAgIHJhbmtPcmkgPC0gYXMubWF0cml4KGdlbmVsaXN0W29yZGVyKC1nZW5lbGlzdFt2aXJ1c05hbWVzW2ldXSksXVt2aXJ1c05hbWVzW2ldXSlcbiAgICAgIHJhbmsgPC0gc3Vic2V0KHJhbmtPcmksIGlzLmZpbml0ZShyYW5rT3JpKSlcbiAgICAgIGdlbmVuYW1lcyA8LSByb3duYW1lcyhyYW5rKVxuICAgICAgXG4gICAgICB0YXJnZXRTZXRzW1tpXV0gPC0gZ2VuZW5hbWVzWzE6dG9wTnVtXVxuICAgICAgXG4gICAgfVxuICB9XG4gIG5hbWVzKHRhcmdldFNldHMpIDwtIHZpcnVzTmFtZXNcbiAgXG4gIGxhYmVsID0gcGFzdGUwKGNlbnRyYWxpdHksIFxcLVxcLCByYW5rTWV0aG9kLCBcXC10b3AtXFwsIHRvcE51bSlcbiAgZGlyLmNyZWF0ZShmaWxlLnBhdGgoXFx2ZW5uXFwsIGxhYmVsKSlcbiAgXG4gIHZlbm5zZXQgPC0gb3ZlckxhcHBlcih0YXJnZXRTZXRzLCB0eXBlPVxcdmVubnNldHNcXClcbiAgdmVubkZpbGUgPC0gcGFzdGUwKFxcdmVubi9cXCwgbGFiZWwsIFxcLnBkZlxcKVxuICBwZGYodmVubkZpbGUsaGVpZ2h0PTYsd2lkdGg9OClcbiAgdmVublBsb3QodmVubnNldCwgbXltYWluID0gcGFzdGUwKFxcSW50ZXJzZWN0aW9ucyBvZiB0YXJnZXQgZ2VuZXM6IFxcLCBsYWJlbCkpXG4gIGRldi5vZmYoKVxuICBcbiAgaW50ZXJMaXN0IDwtIHZlbm5saXN0KHZlbm5zZXQpXG4gIHNhdmUoaW50ZXJMaXN0LCBmaWxlID0gcGFzdGUwKFxcdmVubi9cXCwgbGFiZWwsIFxcLlJEYXRhXFwpKVxuICBcbiAgaW50ZXJOYW1lcyA8LSBuYW1lcyhpbnRlckxpc3QpXG4gIGZvciAoaSBpbiBzZXEobGVuZ3RoKGludGVyTmFtZXMpKSkge1xuICAgIGludGVyTmFtZSA8LSBpbnRlck5hbWVzW2ldXG4gICAgaW50ZXJnZW5lcyA8LSBpbnRlckxpc3RbW2ldXVxuICAgIGlmIChsZW5ndGgoaW50ZXJnZW5lcykgPiAwKSB7XG4gICAgICB3cml0ZS50YWJsZShpbnRlcmdlbmVzLCBmaWxlID0gcGFzdGUwKFxcdmVubi9cXCwgbGFiZWwsIFxcL1xcLCBpbnRlck5hbWUsIFxcLnR4dFxcKSwgc2VwPVxcXFx0XFwsIHF1b3RlPUZBTFNFLCBjb2wubmFtZXMgPSBGQUxTRSwgcm93Lm5hbWVzID0gRkFMU0UpXG4gICAgfVxuICB9XG4gIFxufVxuYGBgXG5gYGAifQ== -->
```r
```r
rankMethod <- \individual\
for (centrality in centralities) {
targetSets <- list()
infile <- sprintf(\./bigTrans-%s-%s-ranks.csv\, centrality, rankMethod)
for (i in seq(length(virusNames))) {
if (file.exists(infile)) {
genelist <- read.csv(infile, header = TRUE, row.names = 1)
colnames(genelist) <- virusNames
rankOri <- as.matrix(genelist[order(-genelist[virusNames[i]]),][virusNames[i]])
rank <- subset(rankOri, is.finite(rankOri))
genenames <- rownames(rank)
targetSets[[i]] <- genenames[1:topNum]
}
}
names(targetSets) <- virusNames
label = paste0(centrality, \-\, rankMethod, \-top-\, topNum)
dir.create(file.path(\venn\, label))
vennset <- overLapper(targetSets, type=\vennsets\)
vennFile <- paste0(\venn/\, label, \.pdf\)
pdf(vennFile,height=6,width=8)
vennPlot(vennset, mymain = paste0(\Intersections of target genes: \, label))
dev.off()
interList <- vennlist(vennset)
save(interList, file = paste0(\venn/\, label, \.RData\))
interNames <- names(interList)
for (i in seq(length(interNames))) {
interName <- interNames[i]
intergenes <- interList[[i]]
if (length(intergenes) > 0) {
write.table(intergenes, file = paste0(\venn/\, label, \/\, interName, \.txt\), sep=\\t\, quote=FALSE, col.names = FALSE, row.names = FALSE)
}
}
}
```